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Abstract 

A  three-dimensional  thermo-fluid-electrochemical  model  is  developed  to  study  the  heat/mass  transport  process  and  performance  of  a  solid  oxide 
fuel  cell  (SOFC).  The  main  objectives  are  to  examine  the  transport  channel  size  effects  and  to  assess  the  potential  of  a  thin-film-SOFC.  A  parametric 
study  was  performed  to  evaluate  the  channel  scale  effects  on  the  temperature,  species  concentration,  local  current  density  and  power  density.  The 
results  demonstrate  that  decreasing  the  height  of  flow  channels  can  lower  the  average  solid  temperature  and  improve  cell  efficiency.  However,  this 
improvement  is  rather  limited  for  the  smallest  channels.  Compared  with  the  conventionally  sized  SOFC,  the  miniaturized  SOFC  with  a  thin-film 
electrolyte  has  the  advantages  of  a  lower  operating  temperature  and  a  better  performance.  Based  on  our  simulation  results,  the  power  density  of  a 
miniaturized  SOFC  could  reach  up  to  5.461  W cm-3.  However,  an  extremely  small  structure  will  lead  to  severe  thermal  stress  induced  by  a  large 
temperature  gradient,  a  cell  with  a  thicker  rib  width  would  have  a  higher  efficiency  and  a  lower  average  temperature.  Numerical  simulation  is 
expected  to  help  optimize  the  design  of  a  solid  oxide  fuel  cell. 

©  2006  Elsevier  B.Y.  All  rights  reserved. 
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1.  Introduction 

Solid  oxide  fuel  cells  are  expected  to  be  widely  applicable  for 
both  small  and  large-scale  power  generation  systems.  The  reason 
is  that  the  solid  oxide  fuel  cell  (SOFC)  is  simple,  highly  efficient, 
tolerant  of  impurities,  and  can  at  least  partially  reform  hydro¬ 
carbon  fuels  internally.  However,  the  application  of  a  SOFC  is 
bounded  by  some  limitations,  which  are  primarily  due  to  its 
high  operating  temperature  (1000-1200  °C),  expensive  materi¬ 
als,  difficulties  in  fabrication,  large  thermal  stresses  within  the 
cells,  etc.  [1,2].  Major  efforts  have  been  placed  on  the  mate¬ 
rials,  designs  of  a  novel  geometry,  and  improvement  in  power 
density  along  with  lowering  the  operating  temperature  down  to 
500-700  °C  [3-10] .  Ramakrishna  et  al.  [4]  developed  a  new  type 
of  SOFC,  which  has  a  thin-wall  geometry.  Because  of  a  better 
longitudinal  reactant  distribution  along  the  thin-wall  geometry, 
they  were  able  to  improve  the  SOFC  performance  and  achieve  a 
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maximum  power  density  of  1.18  W cm-2.  Recently,  Hwang  et 
al.  [6]  have  performed  computational  studies  on  aMOLB  (mono¬ 
block-lay  er-built)-type  SOFC.  In  their  study,  it  is  reported  that 
the  MOLB-type  SOFC  has  a  higher  fuel/oxidant  utilization  ratio 
than  a  planar  SOFC.  Similar  computational  studies  have  also 
been  carried  out  by  Lu  and  Schaefer  for  flat  tube  HPD-SOFC 
[7]. 

In  order  to  minimize  the  overall  losses,  alleviate  the  thermal 
stress  and  improve  the  electrolyte- supported  SOFC  efficiency, 
key  parameters  governing  the  polarizations  and  heat/mass  trans¬ 
fer  efficiency  should  be  investigated.  Basically,  for  a  given  oper¬ 
ating  condition,  there  are  three  kinds  of  mechanisms  that  cause 
the  overall  losses  [1]:  (1)  ohmic  polarization,  which  is  from 
the  ohmic  resistances  of  all  components  and  contact  interfaces; 
(2)  activation  polarizations  at  the  anode  and  cathode,  which  are 
related  to  charge  transfer  processes  and  depend  on  the  nature  of 
electrode-electrolyte  interfaces;  (3)  concentration  polarization, 
which  is  associated  with  the  transport  of  gas  species  and  thus 
depends  on  the  channel  dimensions,  electrode  geometry,  and  the 
nature  of  electrode  microstructure,  i.e.  porosity  and  tortuosity. 
For  an  electrolyte- supported  SOFC,  a  larger  ohmic  polariza- 
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Nomenclature 


A  area  active  area  (m2) 

c  specific  heat  capacity  (J  kg- 1  K- 1 ) 

Ci  molar  concentration  of  component  i  (mol  m-3) 
Dij  binary  mass  diffusivity  (m2  s-1) 

D/  k  effective  Knudsen  diffusion  coefficient  for  com¬ 
ponent  i  (m2  s-1) 

E  Nemst  potential  (V) 

F  Faraday’s  constant,  96486.7  (Cmol-1) 

F\  Forchheimer  coefficient 

A  G  chemical  potential  ( J  mol  - 1 ) 

h  heat  transfer  coefficient  (W  m-2  K- 1 ) 

hmJ  mass  transfer  coefficient  of  component  i  (ms-1) 

H  channel  height  (m) 

A //reform?  A//Shift,  A//eiec  enthalpy  change  of  reactions 
(J  mol-1) 


i 

*oa?  ^oc 


L 

J 

k 


k+ ,  kY  ,k 


K 


Mi 

nQ 

Pi 

Po2 

P 

q " 

^elec 
2oh?  2re 

f 


local  current  denstiy  (Am-2) 

exchange  current  density  of  anode  or  cathode 

(A  m-2) 

channel  length  (m) 

diffusion  mass  flux  vector 

thermal  conductivity  (W  m-1  K-1) 

sh’  ^sh  velocity  constants  for  reforming  or  shift 

reaciton 

permeability  of  the  porous  medium  (m2) 
molar  consumption  or  production  rate  of  com¬ 
ponent  i  at  the  electrode/electrolyte  interface 
(mol  m  2  s-1 ) 

molar  weight  of  component  i  (kg  mol-1) 

number  of  electron 

partial  pressure  of  component  i  (atm) 

partial  pressure  of  oxygen  at  reaction  site  (atm) 

pressure  (Pa) 

heat  flux  at  interface  (W  m-2) 
heat  flux  from  electrochemical  reaction  (W  m-2) 
volumetric  heat  generation  (W  m-3) 
aveage  pore  radius  (m) 


Rax,  Raz,  Rex,  RCz ,  Re,  Rinb,  Rim  resistance  of  electrodes, 
electrolyte  and  interconnectors  (Q) 

Rg  universal  gas  constant,  8.3143  (J  mol-1  K-1) 

Rr,i,  Rsh,i  volumetric  reaction  rates  of  component  i  in 
reforming  or  shift  reaction  (molm-3  s-1) 

Si  net  rate  of  production  for  component  i 
T  temperature  (K) 

V  velocity  vector  (ms-1) 

Xi  molar  fraction  of  component  i 

Yi  mass  fraction  of  component  i 


Greek  letters 

/3  transfer  coefficient 

s  porosity 

act, an?  *7act,ca  activation  overpotential  (V) 
r/con  concentration  overpotential  (V) 

li  viscosity  (m2  s-1) 

Vint  interface  velocity  (ms-1) 


cp  ratio  of  porosity  to  tortuosity 

p  density  (kg  m- 3 ) 

a  electrical  conductivity  ( Q ~ 1  m- 1 ) 

0  electrical  potential  (V) 

Subscripts 
an  anode 

b  bulk 

ca  cathode 

eff  effective 

f  flow 

i  components  in  gas  mixture  (fuel  channel  mixture: 

i  =  CH4,  CO,  CO2,  H2,  H2O;  air  channel  mixture: 
i  =  O2?  N2) 
int  interface 

s  solid 


tion  in  the  electrolyte  is  a  serious  problem,  so  developing  a 
thin-film  electrolyte  to  substantially  reduce  the  ohmic  loss  is  a 
worthwhile  approach  [8-10].  In  general,  with  the  dimension  of 
flow  channel  decreasing  from  conventional  to  microscale  dimen¬ 
sions,  the  temperature/concentration  gradient  decreases  and  the 
heat/mass  transfer  coefficient  increases.  Moreover,  minimized 
sizes  are  important  to  developing  low  cost  and  portable  SOFCs. 
One  therefore  expects  that  concentration  loss  and  operating  tem¬ 
perature  in  smaller  fuel  cells  would  be  effectively  reduced. 

However,  only  a  few  investigations  address  the  issues  of  ther¬ 
mal  stress,  the  heat/mass  transfer  efficiency  and  cell  performance 
for  flow  channels  whose  characteristic  dimension  is  less  than 
1.0  mm  [10-12].  On  the  other  hand,  a  small  scale  may  lead  to 
thermal  management  difficulties  and  poor  thermal  cyclability, 
etc.  Recently,  Li  et  al.  [11]  proposed  a  new  design  for  the  fuel 
cell  gas  distributor  by  maximizing  the  access  area  of  reactants  for 
enhancing  mass  transfer,  and  were  able  to  increase  the  maximum 
output  power  by  40%  for  a  single  cell.  These  potential  advan¬ 
tages  motivate  us  to  investigate  the  scaling  effects,  especially 
of  miniaturized  flow  channels,  on  the  performance  of  a  SOFC. 
Other  geometrical  parameters,  such  as  the  rib  width  and  elec¬ 
trolyte  thickness,  also  have  strong  effects  on  the  performance  of 
a  fuel  cell  and  should  be  optimized. 

Considering  the  complexity  of  the  actual  energy  transport 
process  in  SOFCs,  a  mathematical  model  is  usually  a  desirable 
tool  to  analyze  the  phenomena  inside  a  SOFC  and  predict  its 
performance.  Over  the  past  decades,  various  models  have  been 
proposed  based  on  different  complexity  levels,  ranging  from 
the  molecular-level  via  the  cell  component-level,  and  the  cell- 
level  to  the  system-level  [6,7,13-20].  In  the  present  study,  a 
3D  thermo-fluid-electrochemical  model  was  developed  to  sim¬ 
ulate  the  heat/mass  transport  and  electrochemical  process.  The 
geometry  of  the  interconnector  (or  rib)  and  microstructures  of 
electrodes  were  considered.  Based  on  this  model,  the  effects 
of  channel  height,  thin-film  electrolyte,  and  rib  thickness  on 
the  performance  characteristics  of  a  SOFC  were  quantitatively 
investigated  and  the  following  issues  were  addressed:  (1)  To 
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Fig.  1.  Schematic  of  a  co-flow  planar  SOFC  and  the  computational  domain. 


what  extent  does  the  channel  dimension  affect  the  heat/mass 
transfer?  (2)  As  the  height  of  a  flow  channel  is  decreased, 
whether  or  not  the  power  density  is  improved  and  the  solid  tem¬ 
perature  rise  is  suppressed?  (3)  Whether  or  not  the  thin-film 
electrolyte  improves  the  performance  of  a  SOFC?  (4)  What  is 
the  optimal  thickness  of  the  rib?  (5)  How  do  the  above  variations 
in  geometry  affect  the  thermal  stress? 

2.  Description  of  modeling 

The  stack  of  a  SOFC  is  made  up  of  many  identical  unit  cells 
connected  in  parallel,  most  of  which  function  under  similar  oper¬ 
ation  conditions.  Therefore,  to  reduce  computational  work,  only 
one  unit  cell  is  taken  as  the  computational  domain.  This  typi¬ 
cal  single  cell  is  assumed  to  be  an  interior  and  non-boundary 
unit.  As  a  result,  the  heat  loss  to  the  surroundings  is  not  con¬ 
sidered  and  the  only  heat  removal  mechanism  from  the  cell  is 
through  the  gas  streams.  The  schematic  of  the  computational 
model  is  shown  in  Fig.  1.  Furthermore,  only  one  half  of  the 
gas  channels  are  modeled  by  applying  the  symmetry  condi¬ 
tion.  A  single  cell  consists  of  a  porous  layer  of  anode,  a  porous 
layer  of  cathode,  a  dense  layer  of  electrolyte,  and  interconnec¬ 
tors  (or  ribs).  The  interconnectors  support  the  structure  and  also 
serve  as  the  current  collectors.  When  air  and  fuel  streams  flow 
along  the  X-direction  (co-flow  arrangement),  O2  gas  from  the  air 
stream  diffuses  through  the  cathode  and  accepts  electrons  from 
the  external  circuit  at  the  cathode/electrolyte  interface.  Also,  at 
the  anode/electrolyte  interface  the  oxygen  ions  pass  through  the 
electrolyte  and  react  with  the  H2,  and  the  reaction  product  H2O 
diffuses  back  into  the  fuel  stream.  If  the  fuel  is  natural  gas  (CH4, 
CO2,  CO,  H2,  and  H2O),  the  reforming  and  shift  reactions  take 
place  through  the  anode  layer  as  follows: 

Reforming  :  CH4  +  H2O  3H2  +  CO  +  206kJ/mol 
Water-gas  shift :  CO  +  H2O  H2  +  CO2  —  41  kJ/mol. 


current  are  imposed  at  the  boundaries  of  model  domain.  Walls 
at  the  periphery  of  the  single  cell  are  assumed  to  be  adiabatic. 

2.7.  Thermal- fluid  model 

In  this  model,  the  entire  computational  domain  consists  of 
three  sub-domains:  porous  (electrodes),  solid  (interconnectors 
and  electrolyte)  and  pure  fluid  regions  (channel  flows).  It  has 
been  demonstrated  by  previous  experiments  that  the  fluid  flow  in 
the  porous  media  is  a  low  speed  mixed  gas  flow  with  a  Reynolds 
number  much  less  than  unity.  Therefore,  this  flow  is  treated  as 
laminar  and  incompressible.  It  is  also  assumed  that  the  porous 
medium  is  homogeneous  and  in  local  thermodynamic  equilib¬ 
rium  with  the  fluid.  In  the  porous  region,  the  flow  is  governed  by 
the  Brinkman-Forchheimer  extended  Darcy  model,  which  takes 
into  account  the  effects  of  flow  inertia  as  well  as  friction  induced 
by  macroscopic  shear  stresses.  While  in  the  fluid  region,  the 
Navier- Stokes  equations  describe  the  flow  behavior.  The  con¬ 
servation  equations  of  mass,  momentum,  energy  and  species  for 
fluids  in  both  pure  fluid  region  and  porous  domain  are: 

V  •  (pfV)  =  0  (1) 


— V  •  (VV) 

£ 


S/P 
—  + 

£ 


\M_ 

£ 


V-x 


V  •  (£pfCfVT)  —  V  •  (fceffVr)  +  (2 re  +  <2 oh  (3) 

V  ■  (£PfYiV)  = -V  ■  Ji  +  Si  (4) 

n 

Yr  =  1  (5) 

i= 1 

P(  =  RgT  {Eji/Mi)  (6) 


There  are  several  assumptions  that  should  be  pointed  out 
in  the  present  model:  (1)  transport  and  gradients  in  all  three- 
directions  are  considered;  (2)  energy  transfer  due  to  species 
diffusion  in  channel  flows  and  porous  layers  are  neglected;  (3) 
because  the  single  cell  model  represents  a  repeating  cell  unit 
in  the  center  of  a  larger  stack,  cyclic  boundary  conditions  for 


where  the  parameter  x  is  set  to  unity  for  the  flow  in  a  porous 
medium  (0  <  £  <  1)  and  to  zero  in  regions  without  a  porous  mate- 
rial  (£=  1).  V  and  7/  represent  velocity  vector  and  local  mass 
fraction  of  each  component,  respectively,  s  is  the  porosity  and  K 
is  the  permeability  of  the  porous  medium.  F\ ,  the  Forchheimer 
coefficient  is  calculated  by  F\  =  1. 8/(1 80s5)0-5.  The  effective 
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thermal  conductivity  keff  is  equal  to  sk f  +  (l  —  s)ks.  The  heat 
source  Q0 h  is  due  to  ohmic  loss  and  activation  overpotentials. 
Another  source  term,  gre,  denotes  the  heat  produced  by  reform¬ 
ing  and  shift  reactions,  which  only  takes  place  in  the  anode  layer. 
There  is  no  heat  generation  term  in  pure  fluid  region.  Thus,  once 
the  information  of  the  electrical  field  and  chemical  reactions  is 
obtained,  the  two  heat  source  terms  can  be  exactly  determined. 
Si  is  the  net  rate  of  production  of  component  i  by  chemical  reac¬ 
tion  and  is  left  for  specific  discussion  later.  Ji  is  the  diffusion 
mass  flux  vector  derived  from  the  Stefan-Maxwell  relation  [21] : 


j  +  i 


where  X;  is  the  mole  fraction  of  component  i  and  Ay,eff  is  the 
effective  binary  diffusion  coefficient  of  a  mixture  of  components 
i  and  j  and  equal  to  (pDy.  cp  is  the  ratio  of  porosity  to  tortuosity. 
Ak  is  the  effective  Knudsen  diffusion  coefficient  of  component 
i  and  calculated  according  to  the  kinetic  theory  of  gas: 


rate  due  to  chemical  reaction  in  Eq.  (3)  is  expressed  by: 

Qre  —  AAeformA,CH4  +  A//Shift  7?sh,H2  (16) 

Note  that  Eq.  (16)  is  only  applicable  to  the  anode  layer.  In  Eq.  (3), 
the  second  source  term  includes  the  heat  generations  from  ohmic 
loss  in  all  solid  or  porous  regions  and  activation  overpotentials 
in  electrodes: 

<2oh  =  rR ,  for  solid  or  electrodes  (17a) 

(2  oh  =  rR  -f  r]actI,  for  electrodes  (17b) 

where  I  is  local  current  and  R  is  resistance,  r/act  is  activation 
overpotential.  It  is  assumed  that  the  electrochemical  reaction 
occurs  in  the  vicinity  of  the  anode/electrolyte  interface,  and  so 
the  thermodynamic  heat  generation  from  this  reaction  is: 

q'' iec  =  (A Alec  -  A G)m"2  (18) 

here  AG  is  the  chemical  potential  and  m ^  is  the  molar  con¬ 
sumption  rate  of  hydrogen.  In  the  actual  calculation,  Eq.  (18)  is 
treated  as  a  heat  flux  condition  at  the  anode/electrolyte  interface. 


2.2.  Electrochemical  model 


f  is  the  average  pore  radius.  The  density,  viscosity,  specific  heat, 
and  thermal  conductivity  are  functions  of  both  temperature  and 
mass  fraction  of  gas  mixture.  In  a  SOFC,  the  shift  reaction, 
which  takes  place  at  high  temperatures,  can  be  assumed  to  react 
very  quickly  and  almost  in  equilibrium  [22].  The  reforming 
kinetics  is  assumed  to  be  quasi-homogeneous.  Lehnert  et  al. 
[22]  reported  experimental  results  for  the  forward  velocity  con¬ 
stants  (k+ ,  k+h)  and  the  backward  velocity  constants  (kf ,  k~h )  for 
the  two  reactions  at  three  different  temperatures.  Habeman  and 
Young  [23]  have  developed  empirical  correlations  for  velocity 
constants  with  Arrhenius  curve  fitting  technique  using  exper¬ 
imental  data  of  Lehnert  et  al.  [22]  in  the  temperature  range, 
900  <  T<  1300  K.  The  volumetric  reaction  rates  of  hydrogen  in 
a  shift  reaction  and  methane  in  reforming  reaction  are  calculated 


by: 

Ah,H2  =  k+(RgT)2CcoCu2o  ~  k~ (RgT)2 Cco2Cn2  (9) 

A,ch4  =  (RgT)2Ccn4Cn2 o  -  k~(RgT)4CcoC^2  (10) 

where  C\  (/  =  CH4,  CO,  CO2,  H2,  H2O)  is  the  molar  concen¬ 
tration  of  component  i.  The  net  rate  of  production  for  each 
component,  Su  can  be  stated  as  follows: 

Scu4  =  -Mch4A,ch4  (11) 

Sco  =  A:oA,ch4  -  A:oAh,H2  (17) 

5h2o  =  -Mh2oA,ch4  -  At2oAh,H2  (13) 

Su2  =  3Mh2  A,ch4  +  At2  Ah,H2  (14) 

Sco2  =  Mqq2  Ah,H2  •  (15) 


When  the  molar  numbers  of  methane  and  hydrogen  involved 
in  the  reaction  are  determined,  the  volumetric  heat  generation 


The  electrochemical  model  predicts  the  local  electrical  poten¬ 
tial  and  current  density,  which  is  subject  to  changes  in  geometry, 
local  temperature  and  gas  compositions.  Numerical  simulation 
for  a  three-dimensional  electrical  field  is  directly  analogous  to 
the  calculation  of  heat  transfer  using  a  finite  volume  method. 
A  separate  anode/electrolyte/cathode  unit  is  modeled  and  the 
discretization  of  equivalent  electric  circuit  is  shown  in  Fig.  2. 
Within  an  electrolyte  layer,  the  current  primarily  flows  from  the 
anode  to  the  cathode.  Electrostatic  potential  within  all  compu¬ 
tational  elements  must  satisfy  the  Laplace  equation: 

V  •  (aV0)  =  0  (19) 

where  o  is  the  electrical  conductivity  and  0  is  the  electrical 
potential.  The  calculation  of  an  electric  field  potential  com¬ 
bines  the  following  contributions:  (1)  ohmic  losses  in  all  the 
conducting  materials,  including  the  electrolyte,  the  electrodes, 
and  the  current  collectors,  (2)  the  contact  resistance  at  interfaces 
and  interconnector  resistance,  (3)  activation  overpotentials  r]acu 
which  is  described  by  the  well-known  Butler- Volmer  expression, 
in  the  anode  and  cathode  layers  [16]: 

.  f  /  PneFrjact\ 

1  =  l°  \eXP  )  ~  6XP 

(20) 


where  F  is  the  Faraday’s  constant.  nQ  is  the  number  of  electrons. 
(3  is  the  transfer  coefficient,  z'o  is  the  exchange  current  density. 
These  parameters  are  given  in  Table  1.  If  ^  is  chosen  as  0.5, 
activation  overpotential  for  the  anode  and  the  cathode  can  be 
written  as  follows: 


0  act,  an 


sinh-1  (  — ^ 

F  \2ioJ 


(21) 
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Fig.  2.  Schematic  of  the  equivalent  electrical  circuit. 
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The  concentration  overpotentials  at  electrodes  are  calculated  as: 

^°"  =  ^ln(j|j)  (23) 

where  po2  and  prQ  are  the  partial  pressures  of  oxygen  at  elec- 
trode/channel  interface  and  on  the  reaction  site,  respectively.  The 
distribution  of  local  current  density  must  be  calculated  repeat¬ 
edly  and  the  calculation  must  be  converged  to  satisfy  the  same 
output  voltage  (Vout)  between  the  upper  and  lower  surfaces  of 
the  interconnectors.  Once  the  electrical  potential  at  each  node  is 
known,  the  local  current  density  can  be  calculated  as: 


where  Re  is  the  electrolyte  resistance  and  Aarea  is  the  active  area. 
0ca  and  0an  are  electrical  potentials  at  cathode/air  and  anode/fuel 
interfaces,  respectively. pu2  and  pn2o  are  partial  pressures  of 
hydrogen  and  water  at  interfaces  between  electrodes  and  elec¬ 
trolyte,  respectively.  Eq  =  1.2723  —  2.7645  x  104r  is  a  function 
of  temperature  at  reaction  location  [18].  The  thermo-fluid  model 
and  electrochemical  model  are  coupled  by  the  boundary  condi¬ 
tions  at  interfaces: 


.//  ^  dCo2 


m02  —  ^02,air  0  +  Cq2  '’int.air 


dy 


(26) 


//  9C  h2o 

Wh20  =  —  ^HiO.fuel — — - h  Ch2Q  Vint, fuel 


dy 


(27) 


i(x,  z) 


E  *7actO)  0ca 


area 


(24) 


mH20  —  — mH2 


(28) 


with  Nernst  potential: 


RoT 

E  =  E0  +  In 
IF 


(  PUlPol 
\  ^H2o 


where  m qo  and  m^Q  are  the  molar  production  or  consump¬ 
tion  rate  of  oxygen  and  water  at  the  cathode/electrolyte  and 
(25)  anode/electrolyte  interfaces,  respectively.  Vint?air  and  Vint,fUei  are 
the  interfacial  velocities  induced  by  local  mass  flux  and  defined 
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Table  1 


Model  input  parameters  for  the  baseline  case 


Parameters  and  conditions 

Value 

Fuel  inlet  temperature  (K) 

1100 

Air  inlet  temperature  (K) 

1100 

Inlet  pressure  (Pa) 

1.01  X  105 

Fuel  inlet  velocity  (ms-1) 

0.6 

Air  inlet  velocity  (ms-1) 

2.5 

Cell  length,  L  (m) 

0.05 

Cell  width,  W  (m) 

3.2  x  icr3 

Rib  width,  Wv  (m) 

0.4  x  10“3 

Channel  height,  H  (m) 

1.0  x  10“3 

Channel  width  (m) 

2.4  X  10“3 

Anode  thickness  (m) 

5  x  10“5 

Cathode  thickness  (m) 

5  x  10“5 

Electrolyte  thickness  (m) 

150 x  10“3 

Electric  resistivity  of  electrolyte  (£2  m) 

2.94  x  10-5  x  exp(  10350/7) 

Ionic  resistivity  for  anode  (£2  m) 

2.98  x  10-5  xexp(- 1392/7) 

Ionic  resistivity  for  cathode  (£2  m) 

8.114  x  10"5  xexp(500/7) 

Interconnect  resistivity  (£2  m) 

1  x  10“3 

Contact  resistance  (£lm2) 

VO 

1 

O 

H 

X 

H 

Density  of  electrodes/electrolyte  (kgm-3) 

4400 

Density  of  interconnector  (kg  m-3) 

5700 

Thermal  conductivity  for  anode/cathode  (W  m-1  K-1) 

12 

Thermal  conductivity  for  interconnector  (W  m- 1  K- 1 ) 

11 

Thermal  conductivity  for  electrolyte  (W  m-1  K-1) 

2.7 

Porosity  (%) 

50 

Tortuosity 

3 

Anode  exchange  current  density  (Am-2) 

6300 

Cathode  exchange  current  density  (Am-2) 

3000 

Permeability  of  the  anode  and  cathode  (m2) 

1.7  x  10"10 

Average  pore  radius  (m) 

1  x  10“6 

as: 


Tint,  air  — 


E,-= 


rri: 


z=02,N2  / 


Pair 


Tint,  fuel  — 


£,= 


// 

m- 


/=H2,H2Q,etc.  / 

Pfucl 


(29a) 


(29b) 


3.  Input  parameters  and  numerical  scheme 

The  SIMPLEC  method  is  applied  to  solve  the  discretized 
equations  of  momentum,  energy,  concentration  and  electrical 
potential.  Air  (79%  N2  and  21%  O2)  is  delivered  to  the  air  chan¬ 
nel.  Fuel  (15.22%  CO,  4.52%  C02,  36.8%  H2,  27.3%  H20,  and 
16.16%  CH4)  is  delivered  to  the  fuel  channel.  The  average  cur¬ 
rent  density  is  set  to  4000  A  m-2.  To  study  the  geometry  effects, 
in  each  computational  case  one  geometric  parameter,  such  as 
the  height  of  channels,  rib  width  or  cell  volume,  is  varied  while 
keeping  the  other  input  parameters  fixed. 

Some  input  parameters  for  the  baseline  case  are  extracted 
from  literature  [15,24,25]  and  listed  in  Table  1.  The  physical 
properties  of  common  materials,  which  are  used  for  electrodes 
and  electrolyte,  strongly  depend  on  temperature.  Consequently, 
the  performance  of  a  SOFC  is  significantly  affected  by  the 
temperature  field.  For  instance,  the  resistivity  of  electrolyte  at 
1200  K  could  be  one  or  two  orders  of  magnitude  smaller  than 
that  at  900  K.  This  temperature  effect  is  included  in  the  present 
model.  Since  electrolyte- supported  cells  are  always  designed 


to  operate  at  over  1 100  K  to  avoid  a  larger  ohmic  resistance  in 
the  electrolyte,  the  inlet  fuel  and  air  temperatures  are  chosen  as 
1 100  K  in  the  present  model.  The  local  mass  transfer  coefficient 
hmj  for  component  /,  heat  transfer  coefficient  h ,  and  average 
temperature  T  are  defined  as  follows: 


hm,i  — 


h  = 


171  i 


6"b,/  Qnt,z 
// 


4 


7b  -  7]nt 


(30) 

(31) 


(32) 


where  q'r  is  the  heat  flux  at  the  interface.  CjnU-  and  C^j  are  the 
molar  concentration  at  the  interface  and  for  the  bulk  flow  of 
component  i,  respectively.  A  is  the  surface  area. 


4.  Simulation  results  and  discussion 


4.1.  Baseline  case 


Fig.  3  provides  the  molar  fraction  profiles  for  different  gas 
components  at  an  average  current  density  of  4000  Am-2.  As 
shown  in  Fig.  3(a),  the  oxygen  concentration  distribution  in  the 
z-direction  is  highly  non-uniform.  The  molar  fraction  of  oxygen 
under  the  rib  is  lower  than  those  in  other  areas.  This  is  because 
oxygen  cannot  diffuse  well  underneath  the  interconnector.  On 
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Fig.  3.  Concentration  contours  for  the  baseline  case:  (a)  oxygen  at  the  cathode/electrolyte  interface,  (b)  hydrogen  at  the  anode/electrolyte  interface,  (c)  methane  at 
the  anode/electrolyte  interface,  (d)  temperature  contour  at  z-  0.001  m  (arrow  denotes  flow  direction). 


the  other  hand,  the  rib  does  not  cause  significant  concentration 
gradients  for  H2  in  the  z-direction  since  the  diffusion  of  H2 
is  about  two  to  three  times  faster  than  those  of  other  compo¬ 
nents  as  Fig.  3(b)  indicates.  Similarly,  the  molar  fraction  of  CH4 
decreases  from  the  fuel  inlet  to  outlet  due  to  the  reforming  reac¬ 
tion  (Fig.  3(c)).  The  concentration  gradients  in  the  y-direction, 
which  is  normal  to  the  electrolyte/anode  interface,  produce  some 
concentration  loss,  which  will  be  discussed  later. 

Fig.  3(d)  shows  the  temperature  contour  on  the  z  —  0.001  m 
plane.  Along  the  flow  direction,  the  temperatures  increase  mono- 
tonically  from  the  inlet  to  the  outlet  due  to  the  contribution  of 
Joule  heating  and  chemical  reaction.  For  the  region  near  the  gas 
inlet,  the  temperature  increases  quickly  and  therefore  causes  a 
larger  thermal  stress  than  in  other  areas.  On  the  air  side,  the  tem¬ 
perature  rise  is  not  as  large  as  that  in  the  fuel  side  due  to  high 
mass  fluxes.  Fig.  4(a)  exhibits  the  current  density  contours.  It 
is  found  that  the  local  current  density  profile  has  a  peak  value 
located  roughly  near  the  lower  left  corner.  The  local  current 
density  near  the  symmetric  line  is  relatively  lower  because  of 
the  longest  current  path.  In  the  region  under  the  rib,  the  current 
density  level  is  higher  and  their  profiles  are  very  smooth. 

Fig.  4(b)  shows  the  calculated  overpotential  distributions.  It 
is  found  that  the  ohmic  loss  across  the  electrolyte  and  electrodes 
(symbol  as  Ele. Ohmic)  represents  the  major  source  of  volt¬ 
age  loss,  followed  by  activation  overpotentials,  interconnector 
ohmic  loss  and  concentration  loss.  The  activation  overpotential 
at  the  cathode  and  the  concentration  loss  are  relatively  small. 
Especially,  the  concentration  polarization  is  on  the  order  of  10~4 
(mV)  and  almost  negligible  in  this  case.  In  addition,  further  cal¬ 
culation  indicates  that  this  polarization  is  much  higher  near  the 
rib  due  to  a  weak  diffusion  process.  Clearly,  the  ohmic  heating 
in  an  electrolyte  contributes  most  to  the  increase  in  temperature. 
Therefore,  the  key  elements  to  be  focused  on  for  improvement 
should  be  the  minimization  of  ohmic  polarization  and  sup¬ 
pressing  temperature  rise  without  significantly  decreasing  the 
cell  performance.  In  the  following  sections,  for  the  purpose  of 


improving  power  density,  the  effects  of  the  geometrical  param¬ 
eters  on  the  SOFC  performance  will  be  analyzed. 

4.2.  Effects  of  gas  flow  channel  dimension 

Reducing  the  channel  sizes  results  in  a  shorter  current  path 
and  also  facilitates  higher  heat/mass  transport  rates.  This  will 


Fig.  4.  Performance  for  the  baseline  case:  (a)  current  density  distribution,  (b) 
average  polarizations. 
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Fig.  5.  Effects  of  the  channel  on  the  mass  transfer  rate  (at  z  =  0.001  m  plane):  (a)  H2  mass  transfer  coefficient  effect  on  the  fuel  channel,  (b)  O2  mass  transfer  coefficient 
effect  on  the  air  channel,  (c)  heat  transfer  coefficient  effect  on  the  air  channel. 


lead  to  changes  in  temperature  and  concentration  distributions, 
and  the  power  efficiency  of  a  SOFC. 

Fig.  5(a)-(c)  give  the  calculated  mass/heat  transfer  coeffi¬ 
cients  on  z  =  0.001  m  plane  for  all  cases.  As  shown  in  the  figure, 
the  channels  with  smaller  heights  experience  greater  enhance¬ 
ment  in  the  mass  and  heat  transfer  coefficients.  The  values  of 
coefficients  decrease  sharply  due  to  the  entrance  effect  in  the 
region  near  the  channel  inlets.  Then,  the  downstream  values 
almost  stay  relatively  constant  along  the  flow  direction.  Mass 
accumulation  in  flow  channels  results  in  a  small  increase  in  mass 
transfer  coefficients.  Moreover,  the  increase  in  heat  transfer  rate 
on  the  air  side  will  help  remove  more  heat  and  improve  the 
cooling  effect.  We  may  conclude  that  higher  mass  transfer  rates 
in  smaller  channels  will  significantly  reduce  the  concentration 
polarization  in  channel  flows  and  improve  cell’s  performance. 

The  average  temperature  profiles  for  interconnectors  of  var¬ 
ious  heights  are  shown  in  Fig.  6(a).  When  the  channel  dimen¬ 
sion  is  decreased,  the  average  solid  temperature  level  steadily 
decreases.  However,  the  temperature  slopes  of  smaller  chan¬ 
nels  are  much  steeper  than  those  of  larger  channels.  This  means 
that  although  the  larger  channels  have  higher  temperature  levels, 
the  more  gradual  temperature  increases  will  alleviate  the  risk  of 


severe  thermal  stress  and  thus  help  ensure  the  structural  integrity 
of  the  cell  components.  However,  higher  temperature  levels  put 
a  more  stringing  requirement  for  materials  and  thus  make  fuel 
cells  susceptible  to  thermo-mechanical  failure.  The  average  cur¬ 
rent  density  distributions  for  different  channel  sizes  are  shown 
in  Fig.  6(b).  It  turns  out  that  the  location  of  the  peak  value  for  the 
local  current  density  shifts  towards  the  downstream  direction  of 
the  flow  as  the  channel  dimension  is  decreased.  Moreover,  the 
current  density  profiles  are  flattened  for  larger  channels. 

In  some  of  previous  research  works  [16],  the  interconnectors 
(or  ribs)  were  assumed  to  have  very  high  electrical  conductiv¬ 
ities  so  that  the  solid  resistances  are  neglected.  However,  the 
values  of  ohmic  resistances  in  interconnectors  could  vary  up  to 
one  order  of  magnitude  for  different  materials  [15].  To  deter¬ 
mine  the  effect  of  the  rib  resistance  on  the  cell  performance,  two 
additional  cases,  in  which  rib  resistances  are  modified  to  R  =  0.2 
and  0.02  Q  cm,  were  investigated.  Calculations  show  that  the 
distributions  of  current  density  and  average  solid  temperature 
are  not  sensitive  to  the  change  of  rib  resistance.  However,  the 
terminal  outputs  strongly  depend  on  the  rib  resistances  as  illus¬ 
trated  in  Fig.  6(c).  For  the  cell  with  the  smallest  rib  resistance, 
R  =  0.02  Q  cm,  the  output  increases  slightly  with  decreasing 
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Fig.  6.  Effects  of  the  channel  on  cell  performance:  (a)  average  solid  temperatures,  (b)  average  current,  (c)  output  voltage  comparison. 


channel  dimension.  As  a  contrast,  outputs  are  improved  remark¬ 
ably  with  the  reduction  in  channel  dimensions  for  the  case  with 
higher  rib  resistances.  It  is  interesting  to  note  that  the  cell  with 
the  smallest  channel  dimension  does  not  have  the  maximal  ter¬ 
minal  output.  A  plausible  reason  is  that  on  one  hand,  a  lower 
temperature  increases  the  ohmic  loss  of  solid  part,  which  is  a 
function  of  the  operating  temperature.  While  on  the  other  hand, 
a  shorter  current  path  and  a  lower  concentration  loss  may  par¬ 
tially  or  totally  counteract  this  ohmic  polarization  increase.  In 
addition,  an  extremely  small  structure  will  lead  to  severe  thermal 
stresses  from  a  large  temperature  gradient.  Comprehensively,  the 
cell  with  the  height  of  0.5  mm  shows  the  best  performance  in  the 
present  study. 

4.3.  Evaluation  of  thin-film  electrolyte 

4.3.1.  Temperature  effects  of  the  baseline  case 

As  mentioned  before,  the  operation  temperature  of  a  SOFC 
needs  to  be  maintained  at  a  high  level,  which  is  essential  to  ion 
conductance  in  the  electrolyte  and  high  electrochemical  reaction 
kinetics.  However,  on  the  other  hand,  higher  temperatures  may 
cause  larger  thermal  stresses  and  electrode  sintering  [1].  As  a 
compromise,  the  temperature  of  a  SOFC  is  usually  controlled 
within  a  reasonable  range.  Fig.  7(a)  compares  the  temperature 


profiles  and  ohmic  polarizations  under  two  different  inlet  gas 
temperatures  for  the  baseline  case.  Compared  with  a  96  K  rib 
temperature  rise  in  the  flow  direction  for  the  case  of  a  higher 
inlet  gas  temperature  ( 1 1 00  K) ,  a  larger  temperature  rise  of  1 30  K 
was  found  for  the  case  of  a  lower  inlet  gas  temperature  (1000  K). 
Therefore,  it  is  noted  that  the  temperature  gradient  in  the  ribs 
for  the  lower  inlet  gas  temperature  case  is  30%  higher,  which 
induces  more  thermal  stresses.  The  output  for  cell  with  lower 
inlet  temperature  (1000  K)  decreases  to  0.495  V  from  0.5478  V 
(1 100  K).  In  principle,  lower  flow  temperatures  retard  reforming 
reaction  speed.  In  addition,  it  is  observed  that  the  ohmic  polariza¬ 
tion  is  increased  by  over  50%  due  to  reducing  cell’s  temperature. 
From  the  present  simulation,  it  is  concluded  that  higher  inlet 
temperature  improves  the  efficiency  of  SOFC,  suppresses  the 
temperature  rise  and  thus  avoids  larger  thermal  stress,  however, 
with  the  expenses  such  as  stringent  requirements  for  materials 
and  making  fuel  cell  susceptible  to  thermo-mechanical  failure. 

4.3.2.  Temperature  effects  of  a  thin-film  electrolyte 

One  way  to  regain  a  higher  fuel  cell  performance  at  a  lower 
operating  temperature  is  to  reduce  the  thickness  of  the  elec¬ 
trolyte.  Also,  the  above  results  indicate  that  a  smaller  fuel  cell 
can  greatly  reduce  the  concentration  loss  and  improve  the  fuel 
cell  performance.  Therefore,  a  miniaturized  SOFC  with  a  thin- 
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Fig.  7.  Effects  of  temperature:  (a)  temperature  effect  on  a  conventional  SOFC  (baseline  case),  (b)  the  scaling  effect  on  a  miniaturized-thin- electrolyte  SOFC  with  an 
inlet  temperature  of  850  K,  (c)  efficiency  comparison. 


film  type  electrolyte  is  expected  to  offer  a  better  performance. 
In  the  following  simulations,  the  thickness  of  electrolyte  is  set 
as  10  fju m.  Three  different  cell  volumes  were  investigated,  i.e. 
0.169  cm3  with  a  cross-section  of  1.0  x  3.2  mm2,  0.0448  cm3 
with  a  cross-section  of  0.5  x  1.6  mm2,  and  0.00832  cm3  with  a 
cross-section  of  0.2  x  0.64mm2.  Fig.  7(b)  displays  the  average 
solid  temperatures  and  the  ohmic  polarizations  with  a  reduced 
inlet  temperature  of  850  K.  The  interconnector  temperatures 
decrease  with  decreasing  cell  volumes.  Compared  with  the 
results  for  the  baseline  case  (Fig.  7(a)),  the  total  temperature  rises 
for  the  three  cases  are  still  at  similar  levels,  about  105-135  K.  It 
seems  that  the  thermal  stress  levels  due  to  temperature  gradients 
would  not  be  reduced  for  a  miniaturized  SOFC.  For  the  smallest 
cell,  the  ohmic  polarization  is  higher  than  those  of  the  larger  cells 
due  to  a  lower  operating  temperature.  Fig.  7(c)  demonstrates  a 
miniaturized  SOFC  effectively  improves  the  output  voltage  and 
the  power  density.  For  smallest  fuel  cell,  the  power  density  is 
high,  and  up  to  5.461  W cm-3.  Similar  results  were  reported  in 
recent  work  by  Cha  et  al.  [12].  They  found  that  microchannels 
enhance  fuel  cell  performance  in  terms  of  both  power  density 
and  efficiency.  Experimental  results  further  revealed  that  the 
cell  performance  peaks  for  the  100  jxm  channels  and  gradually 
decreases  with  further  reductions  in  the  cell  size. 


4.4.  Optimization  of  rib  width 

When  the  width  of  the  rib  is  decreased  for  a  fixed  cell  width, 
i.e.  3.2  mm,  the  resistances  in  the  rib  and  the  active  reaction  area 
will  increase.  At  the  same  time,  the  concentration  overpotential 
underneath  the  rib  and  the  contact  resistances  will  decrease.  To 
achieve  the  maximum  cell  performance,  the  rib  width  should 
be  optimized.  Fig.  8(a)  displays  the  effect  of  rib  width  on  the 
temperature  field.  Clearly,  a  wider  rib  reduces  the  rate  of  temper¬ 
ature  increase.  For  example,  total  temperature  rise  is  suppressed 
by  about  30  K  for  the  widest  rib  (0.8  mm)  as  compared  with 
the  thinnest  rib  (0.2  mm).  Thus  it  can  be  concluded  that  a  wider 
rib  can  alleviate  thermal  stresses.  Fig.  8(b)  compares  the  out¬ 
let  molar  fractions  of  various  gas  components  and  the  terminal 
output  for  different  rib  widths.  Having  the  highest  average  solid 
temperature,  the  reforming  and  shift  reactions  are  found  to  be 
most  active  in  the  cell  with  thinnest  rib,  consequently  less  H2 
and  more  CH4  are  consumed  inside  it.  Moreover,  Fig.  8(c)  and 
(d)  indicate  that  for  cells  with  thinner  ribs,  the  concentrations  of 
oxygen  are  more  uniformly  distributed  across  the  whole  porous 
layer  (^-direction),  which  eliminates  inefficiencies  that  would 
normally  be  caused  by  obstructed  “dead  zones”  under  the  ribs. 
Comparing  with  a  thicker  rib,  the  molar  fraction  of  oxygen  near 
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Fig.  8.  Effects  of  rib  width  on  the  cell  performance:  (a)  average  top  interconnector  temperature,  (b)  comparison  of  performance,  (c)  oxygen  at  the  cathode/electrolyte 
interface  with  the  thinnest  rib  (ratio  1:7),  (d)  oxygen  at  cathode/electrolyte  interface  with  the  thickest  rib  (ratio  1:1). 


the  diffusion  layer  surface  is  significantly  higher  for  the  cell  with 
a  thinner  rib.  Similar  tendencies  are  found  for  other  reactants. 
However,  the  decrease  in  concentration  loss  for  a  thinner  rib 
is  not  strong  enough  to  counteract  the  increase  in  ohmic  loss 
and  thus  the  terminal  output  still  decreases.  So  it  is  expected 
that  a  smaller  ratio  of  the  channel  width  to  the  rib  width  would 
serve  to  reduce  the  ohmic  losses  at  the  interface  and  improve 
cell  performance. 

5.  Conclusions 

A  3D  thermo-fluid/electrochemical  model  was  developed 
to  evaluate  the  heat/mass  transfer  and  electrochemical  perfor¬ 
mance  for  an  electrolyte- supported  SOFC.  Based  on  this  model, 
detailed  investigations  were  carried  out  to  study  the  effects  of 
operating  temperature  and  geometry  on  the  distributions  of  tem¬ 
perature,  local  concentration,  local  current  density  and  output 
voltage.  The  main  conclusions  reached  are  listed  below: 

(1)  When  the  height  of  flow  channels  is  decreased,  the  average 
solid  temperature  is  reduced  accordingly  and  the  cell  effi¬ 
ciency  is  improved  due  to  both  a  higher  heat/mass  transfer 
coefficient  between  the  channel  wall  and  flow  stream  and  a 
shorter  current  path.  However,  a  larger  fuel  cell  has  a  more 
gradual  temperature  rise,  which  alleviates  the  risk  of  severe 
thermal  stress  and  helps  ensure  the  structural  integrity  of 
the  cell  components.  If  the  channel  height  is  shorter  than 
0.5  mm,  then  the  cell  efficiency  will  not  be  improved  any 
further. 


(2)  For  a  conventional  SOFC,  if  the  gas  inlet  temperature  is 
higher,  then  higher  cell  efficiency  will  result,  which  causes 
an  increased  cell  terminal  output,  a  lower  solid  temper¬ 
ature  rise  across  the  channel  and  smaller  thermal  stress 
accordingly.  However,  highly  elevated  solid  temperatures 
usually  make  significant  demands  on  materials.  In  compar¬ 
ison,  the  miniaturized  SOFC  with  a  thin-film  electrolyte 
has  the  advantages  of  a  lower  operation  temperature  and 
an  excellent  performance.  Conservatively,  we  expect  that 
the  power  density  of  a  miniaturized  SOFC  with  a  thin-film 
electrolyte  will  have  the  potential  to  exceed  lOWcm-3. 
Of  course,  the  average  solid  temperature  may  be  relatively 
high  and  therefore,  thermal  management  will  be  a  challenge. 
However,  it  is  noted  that  if  the  radiation  loss  to  the  environ¬ 
ment  is  considered,  the  overall  energy  loss  may  be  severe 
and  the  performance  will  be  lowered  because  of  a  larger 
surface  to  volume  ratio.  In  addition,  the  thermal-mechanical 
reliability  of  thin  films  should  be  considered  in  actual  appli¬ 
cations  [8-10]. 

(3)  The  smaller  ratio  of  the  channel  width  to  the  rib  width  con¬ 
tributes  to  a  reduction  of  the  ohmic  loss  at  the  interface  and 
thus  improves  the  cell  performance. 

Overall,  the  present  model  has  demonstrated  that  there  is  a 
trade-off  among  competing  factors  influencing  the  performance 
of  a  SOFC,  which  was  not  fully  addressed  in  previous  work. 
It  should  also  be  noted  that  reducing  the  cell  volume  and  the 
electrolyte  thickness  and  optimizing  the  rib  width  are  only 
some  of  the  factors  influencing  the  performance  of  a  SOFC. 
Therefore,  finding  better  approaches  or  other  critical  parameters 
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to  optimize  the  design  of  a  SOFC  will  be  a  subject  for  future 
work. 
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